%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                                                                         %
%              Population Growth and Firm Dynamics                        %
%                Authors: Michael Peters and Conor Walsh                  %
%                                                                         %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


% ---------------------------- Description -------------------------------- 
    
    % This is the main file of the code. By running this script, the user
    % performs the following tasks:
        
        % (a) Calibration: the user finds the set of parameters that match
        % the model's moments to the moments that are being targeted on the
        % data.
        
            % By default, Compute_Calibration = 0 (Master, Line 90).
            % This means that the parameters are loaded directly from
            % paramguess.mat. If Compute_Calibration = 1, then the code
            % will search for the calibration parameters.

        % (b) Model: the user uses a parameter vector to compute the main
        % objects of the model.

        % (c) Counterfactual: the user computes the main objects of the
        % model under a counterfactual populational growth rate.

        % (d) Transitional Dynamics: with the model's equilibrium
        % conditions, the user derives the dynamic response of the model
        % from the original BGP to the counterfactual BGP.

            % By default, the variable Compute_TD = 0 (Master, Line 147).
            % This means that the output is loaded directly from
            % InfoMatrix.mat. If the user wants to compute the Transitional
            % Dynamics, Compute_TD = 1 needs to be set.

% -------------------------- Important Notice -----------------------------  

    % Whenever a comment refers to an Equation, it will be according to the
    % following standard:

        % Equation (M-X) - Equation X on the Main Article;
        % Equation (A-X) - Equation X on the Appendix;
        % Equation (SM-X) - Equation X on the Supplementary Material;
        % Equation (RP-X) - Equation X on the Replication Package;

    % This same terminology will be applied to refer to Sections and others.

%--------------------------------------------------------------------------
%                              Parameters                                 %
%--------------------------------------------------------------------------

    % Clears the Workspace
    clear

    % Set invariant parameters of the model
    Params

%--------------------------------------------------------------------------
%                            Target Moments                               %
%--------------------------------------------------------------------------

    % Aggregate Productivity Growth Rate
    p.g_Moment = 0.0184;

    % Markup Growth by Age 10 (from RevLBD)
    p.MarkupGrowth_Moment = 0.102; 

    % Sales Growth by Age 10 (from RevLBD)
    p.SalesGrowth_Moment = 0.58;    

    % Average Firm size in Employees per Firm (from BDS)
    p.AverageFirmSize_Moment = 19.816; 

    % Aggregate Entry Rate in 1980 (from BDS)
    p.EntryRate_Moment = 0.116;             

    % Average Cost-Weighted Markup in 1980
    p.Markup_Moment = 1.20;     

    % Employment share of firms with +10,000 employees
    p.ShareLarge_Moment = 0.262;

%--------------------------------------------------------------------------
%                              Calibration                                %
%--------------------------------------------------------------------------

    % Switch to determine whether Calibration computed (1) or loaded (0)
    Compute_Calibration = 1;    

    % Calibrate the model according to the "Calibration" script
    Calibration;

%--------------------------------------------------------------------------
%                                Model                                    %
%--------------------------------------------------------------------------

    % This scripts computes the model for the calibrated parameters
    Model;

%--------------------------------------------------------------------------
%                            Counterfactual                               %
%--------------------------------------------------------------------------

    % ----------------- Population Growth Trajectory ----------------------

    % Import file in which we have the sequence for population growth
    etaseq = csvread('LF_Growth_QuantAnalysis.csv',1,0);

    % We make sure the series starts on 1981 (after year of calibration)
    etaseq = [etaseq(5:end,1)/100, etaseq(5:end,2)];

    % We add our base population rate to the year 1980
    etaseq = [p.eta, 1980; etaseq];

    % We interpolate the series linearly to our time unit
    eta_vector = interp1(etaseq(:,2), etaseq(:,1), p.TimeGrid+1980, 'linear', 'extrap');

    % We smooth the series to avoid very sharp corners
    [eta_vector, ~] = hpfilter(eta_vector, Smoothing = 300);
    eta_vector = eta_vector';

    % Once filtered, we make sure the path starts at the calibrated BGP
    eta_vector(1) = p.eta;

    % -------------------- Model under new BGP ----------------------------

    % Determine the population growth rate under the new BGP
    p.eta = eta_vector(end);

    % Compute counterfactual model under new BGP
    [g_counter, LFN_counter, AvFirmSize_counter, CostMarkup_A_counter, ...
        CostMarkup_counter, EntryRate_counter, Lam_counter, M_counter, ...
        prod_share_counter, LPNwMLam_counter, k_func_lam_counter, k_func_mu_counter, ...
        ExitRate_counter, xProds_counter, ProductDist_counter, Firms_10y_counter, ...
        FirmAge_dens_counter, Exit_age_counter, Esales_counter] = Counterfactual(p);

    % Restore the original population growth rate
    p.eta = eta_vector(1);

%--------------------------------------------------------------------------
%                        Transitional Dynamics                            %
%--------------------------------------------------------------------------

    % Switch to determine whether TD is being computed (1) or loaded (0)
    Compute_TD = 1;

    % Compute the transitional Dynamics
    Transitional_Dynamics;  

    % Computes the main firm moments along the transition path
    Firm_Moments; 


%--------------------------------------------------------------------------
%                               Figures                                   %
%--------------------------------------------------------------------------

    % Plots the main figures displayed on the article
    Figures;

    % Compute Welfare Consequence
    Welfare;

%--------------------------------------------------------------------------
%                     Save the Data For Plotting                          %
%--------------------------------------------------------------------------

% Save the data as CSV files
    SaveFigures;

    